Critical behavior and universality in Levy spin glasses 
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Using large-scale Monte Carlo simulations that combine parallel tempering with specialized cluster 
updates, we show that Ising spin glasses with Levy-distributed interactions share the same univer- 
sality class as Ising spin glasses with Gaussian or bimodal-distributed interactions. Corrections to 
scaling are large for Levy spin glasses. In order to overcome these and show that the critical expo- 
nents agree with the bimodal and Gaussian case, we perform an extended scaling of the two-point 
finite size correlation length and the spin-glass susceptibility. Furthermore, we compute the critical 
temperature and compare its dependence on the disorder distribution width to recent analytical 
predictions [J. Stat. Mech. (2008) P04006]. 

PACS numbers: 75.50.Lk, 75.40.Mg, 05.50.+q 
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I. INTRODUCTION 

Although universality has been established for many 
systems without disorder and frustration, there are still 
skeptics that question this cornerstone of the theory of 
statistical mechanics when applied to disordered spin sys- 
tems with frustrated interactions. According to univer- 
sality, the values of quantities such as critical exponents, 
do not depend on microscopic details of the model, but 
only on e.g., the space dimension and the symmetry of the 
order parameter. Arguments based on high-temperature 
series expansions^ support universality and there is no 
a priori reason why systems with both disorder and frus- 
tration, such as spin glasses^ might not show universal 
features. However, numerical studies are difficult^— and 
suffer from strong corrections to scaling. Therefore, there 
is still debat o n i 16 i 18 i 21 i 29 for some model systems if the 
shape of the disorder distribution can influence the uni- 
versality class of the system. 

Although it is now well established that universal- 
ity is not violated for nearest-neighbor spin glasses 
with compact disorder distributions (e.g., Gaussian or 
bimodal) t 22 i 26 i 27 some studies suggest that this might not 
be the case when the disorder distributions are broad^ 
If the spin interactions are drawn from a Gaussian or 
bimodal distribution the probability to have extremely 
large interactions is very small. It is, however, unclear if 
strong couplings between the spins change the universal- 
ity class of the system. Selecting the interactions between 
the spins from a Levy distribution allows one to contin- 
uously tune the probability to have very strong bonds in 
the system. In particular, for a < 2 (see below for de- 
tails) the Levy distribution has broad tails and thus the 
probability to have a strong bond between two spins is 
large, especially in the limit a — > 1. 

Using large-scale Monte Carlo simulations that com- 
bine parallel tempering with specialized cluster moves^ 
as well as extended scaling techniques,— our results show 
that Levy spin glasses do obey universality for the sys- 
tem sizes studied. Our estimates of the critical ex- 



ponents agree within error bars with the best known 
estimate d 26 ' 27 for Gaussian and bimodal disorder. Fur- 
thermore, we probe recent analytical predictions^ made 
for the critical temperature of Levy spin glasses as a func- 
tion of the disorder distribution width. 

The paper is structured as follows: In Sec. |H] we in- 
troduce the model studied, as well as the measured ob- 
scrvablcs. Section Hill outlines the special (cluster) algo- 
rithm used to treat strong interactions in the Levy spin 
glass, the finite-size scaling analysis, and how we estimate 
the critical temperature, followed by results presented in 
Sec. IIV1 as well as concluding remarks. 




FIG. 1: (Color online) Levy distribution V(J) for the differ- 
ent values of the shape parameter a and c = 1/ \/2, as studied 
here. In particular, for a = 2 a Gaussian distribution is re- 
covered. For 1 ^ a < 2 the distribution is fat-tailed, as can 
be seen in the linear- log plot. 
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II. MODEL AND OBSERVABLES 

We study the critical behavior of the Edwards- 
Anderson Ising spin glass^i with Levy-distributed inter- 
actions, i.e., 



ijSiSj, 



(1) 



where the sites i lie on a three-dimensional cubic lattice 
of size N = L 3 , L the linear dimension, and the spins 
Si can take the values ±1. Periodic boundary conditions 
are used to reduce corrections to scaling. The sum is 
over nearest neighbors and the interactions are inde- 
pendent random variables taken from a Levy distribution 
with zero mean and c = l/\/2 defined through the char- 
acteristic function 4>{t) as 



1 f + °° 
r(J) = ^ / ' 

27r J-oo 



dt<p{t)t 
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(2) 

The parameter a influences the shape of the distribution 
and, in particular, the width of the tails. For a = 2 
Eq. ^ reduces to a Gaussian with variance er = v2c. 
When 1 ^ a < 2 the tail of the distribution decays as 
a power law, as seen in Fig. [T] In this case, exchange 
interactions with very large values can occur, albeit with 
small probability. However, such strong interactions form 
"dimcrs" of spins that cannot be flipped with standard 
Monte Carlo methods^ For decreasing a the probabil- 
ity to have dimers grows, as well as the average of the 
maximum exchange interaction value. 

To test universality, at least two independent critical 
exponents need to be computed. Therefore, in the simu- 
lations, we measure the following quantities. 

The spin overlap q is defined as 



(1) 5 (2) 



(3) 



where (1) and (2) are two copies of the system with iden- 
tical disorder. Using q we define the Binder ratio^ g 
via 



3- 




G[lV"(t-T c )], (4) 



where (• • • ) represents a thermal average and [• • • ] av an 
average over the disorder. The Binder ratio is a dimen- 
sionless function G, i.e., data for different system sizes 
cross at a putative transition temperature T c . A finite- 
size scaling analysis of the universal functional G allows 
one to determine the critical exponent v for the correla- 
tion length. 

The spin-glass susceptibility %SG is defined as 



Xsg = N [(q 2 )] av ~ 1?-*C\I}I\T - T c )} 



(5) 



A finite-size scaling analysis of the susceptibility thus per- 
mits the calculation of the critical exponent 77. However, 



a simple scaling analysis of the spin-glass susceptibility 
suffers from strong corrections to scaling^ and therefore 
an extended scaling^ is performed below where the scal- 
ing function incorporates corrections derived from the 
resummation of a high-temperature series expansion. 

Finally, we measure the two-point finite-size correla- 
tion lcngth.i^ To do so we introduce the wavc-vector- 
dependent spin-glass susceptibility 

XSG(k) = ^[(S iSj )\ v e^-^. (6) 

The two-point finite-size correlation length is then 
given by 



1 



Xsg(O) 



2sin(fe min /2) y xsG(kmin) 
where k m ; n = (27r/L,0,0). It scales as 



- 1, 



L 



(7) 



(8) 



i.e., whenever T = T C data for different system sizes cross 
at one point, up to corrections to scaling £L 



III. NUMERICAL DETAILS 

To test for universal behavior a detailed numerical 
study needs to be performed where one has to ensure 
that the data are in thermal equilibrium. For this pur- 
pose we use a special cluster algorithm that ensures that 
spin dimcrs flip in reasonable simulation times. In addi- 
tion, we describe the data analysis used. 



A. Algorithm 

The simulations are done using the parallel tempering 
Monte Carlo method^ combined with a special cluster 
flip algorithm^ that ensures ergodic behavior even in 
the presence of excessively strong exchange interactions 
between few spins. 

Because the Levy distribution has power-law decaying 
tails, for certain values of the parameter a the exchange 
interactions Jy can be very large. If two spins have a 
strong interaction they will be virtually "frozen" under 
single-spin-flip dynamics. To avoid extremely long equi- 
libration times, at the beginning of each simulation dif- 
ferent sets of clusters C n are generated^ The generation 
of the sets C n is done the following way: 



1. Set jo 



T max /4, where T max is the maximal 



temperature from the simulated temperature set. 

2. The clusters in set C n consist of spins connected by 
bonds that satisfy |Jy| > J™ in . 
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3. The cluster set C„ is stored if C n 7^ C„_i (or if 
n = 0). 

4. J" in is iteratively incremented by one (J™^ 1 = 
J™ in + 1). The procedure is repeated initiating from 
step two until C n consists only of clusters of size 2. 
During the procedure all clusters are stored. 

One Monte Carlo sweep consists of the following proce- 
dure: Each spin of the system is picked once. After hav- 
ing picked the spin, a single-spin flip is performed with 
probability p = 0.75 (empirically we find that for p ~ 0.75 
equilibration is fastest), otherwise a cluster move is done. 
In particular: 

□ The single spin flip is done with the Metropolis 
probability min{l, exp(— AE/T)}, where AE is the 
energy difference between the current configuration 
and the configuration with the spin flipped. 

□ The cluster flip algorithm works as follows: One 
cluster from all sets is randomly (uniformly) 
picked and flipped with the Metropolis probabil- 
ity min{l,cxp(-A£ , /T)}, where AE is the differ- 
ence between the energy of the actual configuration 
and the configuration with the cluster flipped. The 
cluster flip is independent of the orientation of the 
spins in the cluster, i.e., the clusters contain only 
spin indices, such that each spin in the cluster can 
change the direction by other update steps. 

Note that typical cluster sizes range from 2 -20 spins. 

Because the equilibration test for Gaussian disorder- 
does not work when the disorder is Levy distributed, the 
equilibration is monitored by logarithmic binning. All 
measured observables (and their higher moments) are 
recorded as a function of simulation time. Once the last 
four bins agree within error bars the system is deemed 
to be in thermal equilibrium. If this test is not passed, 
the simulation time is increased by a factor of 2 until 
this is the case. Simulation parameters are summarized 
in Table HI 

B. Finite-size scaling analysis 

To gain insights on the strength of the correc- 
tions to scaling, we can compare two dimensionless 
quantities ; 22 ' 26 the correlation length £l/L and the 
Binder parameter g. By plotting <7[£z,(T, L)/L] there are 
no nonuniversal metric factors. Therefore, data for all 
system sizes simulated and a given parameter a should 
all collapse onto a universal function if there are no cor- 
rections to scaling. Data for a = 1.25 are shown in Figure 
[2] and illustrate that corrections are large for L < 8. 

Furthermore, if two different models share the same 
critical exponent v, because no nonuniversal factors when 
plotting g[^L(T,L)/L] are present, all data should col- 
lapse onto a universal curve. In Fig.[2]we also show data 
for Gaussian disorder (a = 2) for a large system size 



TABLE I: Parameters of the simulations for different a val- 
ues. iV sa is the number of samples, iV sw is the total number of 
Monte Carlo sweeps used for equilibration (the same amount 
is used for measurement), T m i n is the lowest temperature sim- 
ulated, T max is the highest temperature simulated, and Nt is 
the number of temperatures used in the parallel tempering 



method for each system size L. 


a L 7V sa 7V SW 


T ■ 

± mm 


^max 





1.00 
1.00 
1.00 
1.00 
1.00 
1.00 

1.25 
1.25 
1.25 
1.25 
1.25 

1.50 
1.50 
1.50 
1.50 
1.50 

1.75 
1.75 
1.75 
1.75 
1.75 



(L = 24) M Data for a = 1.25 and L > 8 agree with the 
Gaussian case, thus illustrating that for a conventional 
scaling analysis only the largest system sizes should be 
included. 

0.9 , , , , , 1 




ih{T,L)/L 

FIG. 2: (Color online) Binder ratio g as a function of the 
finite-size correlation length £l/L for several system sizes and 
for a Levy parameter a = 1.25. Strong corrections to scaling 
are visible. The data for the largest system sizes simulated 
agree with the Gaussian case (a = 2, L = 24, from Ref . l26h . 

We have attempted different scaling approaches , 22 ' 26 
as well as the inclusion of corrections to scaling^ How- 
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ever, large system sizes arc difficult to simulate for Levy 
spin glasses and therefore we use the extended scaling 
technique 2 ^ that allows us to include smaller system sizes 
in the scaling analysis. 

Within the extended scaling framework 2 ^ the standard 
scaling expression for the correlation length, Eq. ([8]), is 
replaced by 

^ ~ X[(LT) l l v \l - {T/T c f\}. (9) 

The aforementioned expression is derived by including a 
resummation of a high-temperature series expansion and 
therefore includes effects of corrections to scaling. Sim- 
ilarly, the scaling relation for the susceptibility, Eq. 
is replaced by 

Xsg(L,T) ~ (LT) 2 -"C[(LT)V-|1 - (T/T c ) 2 \}. (10) 

We assume that the scaling function in Eq. can be 
approximated by a third-order polynomial for tempera- 
tures larger than T c , i.e., X(x) = a + bx + cx 2 + dx 3 where 
x = (LT) 1 / l '\l - (T/T c ) 2 \, T > T c and perform a fit to 
the six parameters a, b, c, d, T c and v. A similar ap- 
proach is used for the spin-glass susceptibility [Eq. (fT0|)], 
for which there is a seventh parameter, the critical ex- 
ponent r\. Note that the extended scaling scheme only 
works for temperatures T > T c , Therefore, we first per- 
form a rough estimate of T c using conventional scaling 
methods. The nonlinear fit is performed with the statis- 
tics package R^I including system sizes L ^ 6. Error 
bars are determined using a bootstrap analysis. 

To compute the error bars we apply the following pro- 
cedure: For each system size L and iV sa disorder realiza- 
tions, a randomly selected bootstrap sample of _/V sa disor- 
der realizations is generated. With this random sample, 
an estimate of the different observables is computed for 
each temperature. We repeat this procedure A^oot = 500 
times for each lattice size and then assemble A^oot com- 
plete data sets (each having results for every size) by 
combining the i-th bootstrap sample for each size for 
i = 1, . . . , -/Vboot- The finite-size scaling fit described 
above is then carried out on each of these iVboot sets, 
thus obtaining -/Vboot estimates of the fit parameters. Be- 
cause the bootstrap sampling is done with respect to the 
disorder realizations which are statistically independent, 
we can use a conventional bootstrap analysis to estimate 
statistical error bars on the fit parameters. These are 
comparable to the standard deviation among the -/Vboot 
bootstrap estimates^ 

From the aforementioned finite-size scaling we can also 
extract the critical temperature T c (a) to compare to an- 
alytical predictions. We use the critical temperature es- 
timated using the extended scaling method for the corre- 
lation length [Eq. (J9))] because corrections to scaling are 
smaller than for the spin-glass susceptibility. Further- 
more, the bootstrap analysis requires one parameter less 
leading to smaller statistical errors. 



TABLE II: Summary of estimates of the critical parameters. 
T| and iA are the critical parameters estimated from an ex- 
tended scaling analysis of the two-point correlation length, 
whereas r] x has been computed from a finite-size scaling anal- 
ysis of the spin-glass susceptibility with T c , v and 77 as free 
parameters. t\ x \t c .v is computed from a finite-size scaling 
analysis of the susceptibility with T c = and v = iA fixed 
and only 77 as a parameter. 
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r] x 


V x \t c ,v 


1.00 


1.467(31) 


2.42(17) 


-0.346(220) 


-0.438(26) 


1.25 


1.209(28) 


2.49(17) 


-0.411(209) 


-0.412(20) 


1.50 


1.094(21) 


2.34(15) 


-0.344(274) 


-0.414(20) 


1.75 


0.996(20) 


2.61(19) 


-0.274(224) 


-0.413(17) 



IV. RESULTS 

Corrections to scaling for small systems of Levy spin 
glasses with 1 ^ a < 2 are large (see Fig. [5]). We attempt 
to scale the data using the extended scaling scheme, as 
shown in Fig. [3] (center and right columns). The left col- 
umn shows the finite-size correlation length for different 
values of the parameter a. In all cases the data cross at a 
transition temperature that decreases with increasing a. 
In the center panels of Fig. |3]we show an extended scaling 
of the two-point finite size correlation length according 
to Eq. ([9]) with the critical exponents v and T c as param- 
eters. The right column of Fig. [3] shows an extended scal- 
ing of the spin-glass susceptibility according to Eq. (fT0|) 
with 77, v and T c as free parameters. The scaling of the 
data works well and, in particular, the estimated critical 
exponents agree with the bimodal values^ Our best esti- 
mates are summarized in Table HI] Furthermore, in Fig. [31 
we compare our estimates for 77 and v to the bimodal es- 
timates [77 = -0.375(10) and v = 2.45(15)]^! The data 
therefore suggest that all studied Levy spin glasses share 
the same universality class. 

To further strengthen our results for the finite-size cor- 
relation length, in Fig. [5] we show g[^L(L,T)/L] for the 
largest system size studied and different a, as well as 
data for Gaussian disorder^ The data collapse cleanly 
onto a universal curve without any scaling parameters 
providing further evidence for universal behavior. The 
inset of Fig.[5]shows £l/L (T = T c ) for different values of 
the exponent a. For all cases the data agree within error 
bars with the best estimate for bimodal disorder— hence 
strengthening our claim for universal behavior. 

Finally, we show in Fig. [5] estimates for the critical 
temperature T c as a function of the exponent ct i 30 ' 38 The 
horizontal blue line represents the Gaussian limiti 26 ' 27 
The red curve represents T c (a) for a mean-field spin- 
glass model on a diluted graph with fixed connectivity 
k + 1 = 6. The critical temperature is determined from 
the following equation 

1 = k J dJV( J) tanh 2 ( J/T c ), (11) 

where V{J) is given by Eq. There is qualitative 
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FIG. 3: (Color online) Left: Two-point finite-size correlation length £,l/L as a function of the temperature T for different Levy 
parameters a. The data cross, thus signaling the presence of a transition. Center: Extended scaling of the two-point finite-size 
correlation length for different a. Right: Extended scaling of the spin-glass susceptibility for different a. The data scale very 
well and the critical exponents extracted from the scaling agree within error bars, thus suggesting that all systems share the 
same universality class. See Table HT1 for the optimal scaling parameters. 
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FIG. 4: (Color online) Critical exponents r\ and ^ as a func- 
tion of the Levy parameter a. The shaded areas correspond 
to the estimates for Bimodal disorder from Ref. I27I . The es- 
timate for the critical exponent t/* comes from an extended 
finite-size scaling analysis of the two-point correlation length. 
The estimates for the critical exponent r\ are from two inde- 
pendent analyses of the spin-glass susceptibility. rj x is com- 
puted from an extended finite-size scaling analysis where 77, 
v, and T c are parameters, whereas t) x \ v ,t c is computed by fix- 
ing v = jA and T c = T§ from the analysis of the two-point 
correlation length. For all values of a studied, the exponents 
1/* and r) x are in good agreement with the best-known esti- 
mates for the bimodal case. However, the estimate for rf-\ v ,T a 
consistently lies above the best estimate for r\ possibly due to 
strong corrections to scaling that we cannot account for, as 
well as systematic errors from the determination of iA— 



agreement in the trend of the data. At first sight, there is 
a disagreement to the behavior obtained for the infinite- 
range model studied in Rcfs. andl3fl because the de- 
pendence on a is reversed in that case. However, the 
large connectivity limit of Eq. amounts to the ex- 
pression for the critical temperature stated in Refs. 
and [30I For the infinite-range model an a-dependent 
rescaling of the couplings (Jy — > J^N^ 1 ' 01 ) is neces- 
sary to obtain a non-trivial thermodynamic limit which 
changes the energy scale in an a-depcndcnt way. 



h 

1-5 



^.68 

II 0.66 
^"0.64 








-i 


! 

[ 


] 

[ 


n: 

1 " 


^ 0.62 








0.6 0.8 1 

Zl(T,L)/L 

FIG. 5: (Color online) Binder ratio g as a function of the 
two-point finite-size correlation length £l/L for a — 1.00, 
1.25, 1.50 and 1.75 for system sizes L = 10 and 12 (14, for 
a = 1.00), as well as 16 and 24 for the Gaussian (a = 2.0) 
case.— The line is a guide to the eye. All data collapse 
onto a universal curve, thus providing further evidence for 
universality. The inset shows £l / X (T = T c ) — also a univer- 
sal quantity — as a function of a. For all values of a stud- 
ied the data agree within error bars. The horizontal shaded 
area corresponds to the best estimate for bimodal disorder 
£l/L (T = T c ) = 0.645(15).-^ 




V. CONCLUSIONS 

We have studied the critical behavior of a three- 
dimensional Ising spin glass with Levy-distributed inter- 
actions to test universality. An extended scaling analysis 
of the correlation length and spin-glass susceptibility sug- 
gests that for all values of a the Levy spin glass obeys 
universality. Previous claims that universality might be 
destroyed when a — > 1 possibly stem from the fact that 
the simulations did not take into account the effects of 
the strong interactions between some spins, i.e., render- 
ing the simulations noncrgodic. Further support for uni- 



FIG. 6: (Color online) Critical temperature T c computed 
from a finite-size scaling of the two-point correlation length, 
Eq. (J7J, using the extended scaling technique, Eq. The 
continuous curve (labeled with VB) is the critical temperature 
T c for a mean-field spin glass on a diluted graph with fixed 
connectivity k + 1 = 6. The horizontal thick line represents 
the critical temperature T c for the Gaussian spin glass.— 



versal behavior is given by the plot of g against £l/X, 
Fig. where data for all the models studied collapse onto 
a single universal curve. We do find strong corrections 
to scaling and therefore studies with larger system sizes 
and a clear understanding of scaling corrections would 
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be desirable. However, we find no clear evidence for the 
lack of universality. 
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